非線形回帰分析 - 18

指数関数のパラメータの推定

a2,をシフト

\(\Large \displaystyle y_i = a_0 \ Exp (- a_1 x_i) + a_2 \)

 

a0       10.5226 10.4345 10.3454 10.2552 10.1640 10.0718 9.9786 9.8845 9.7894
a1       0.4462 0.4577 0.4696 0.4819 0.4946 0.5077 0.5213 0.5355 0.5501
a2       0.4839 0.5839 0.6839 0.7839 0.8839 0.9839 1.0839 1.1839 1.2839
δ       -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4
                         
i x y   \( \hat{y} \)
1 0 10   11.0065 11.0185 11.0293 11.0391 11.0479 11.0557 11.0625 11.0684 11.0733
2 2 4   4.7943 4.7612 4.7284 4.6960 4.6639 4.6324 4.6015 4.5712 4.5416
3 3 2   3.2427 3.2270 3.2128 3.2001 3.1891 3.1799 3.1724 3.1668 3.1632
4 4 1   2.2496 2.2563 2.2651 2.2763 2.2897 2.3056 2.3239 2.3447 2.3680
5 6 0.5   1.2072 1.2534 1.3021 1.3532 1.4067 1.4627 1.5210 1.5817 1.6447
6 9 0.1   0.6735 0.7535 0.8350 0.9181 1.0025 1.0883 1.1754 1.2637 1.3532
         
S (\(y_i - \hat{y} \)の平方和)       0.4311 0.3554 0.2996 0.2650

0.2531
(Se)

0.2655 0.3035 0.3687 0.4626
dS (Seとの差分)       0.1780 0.1023 0.0464 0.0118 0 0.0123 0.0503 0.1155 0.2094
                         

 

・残差平方和

推定値からの残差

\(\Large \displaystyle Se = \sum_{i=1}^{n} \left[ y_i -\hat{a_0} \ Exp(- \hat{a_1} \ x_i) - \hat{a_2} \right]^2 \)

a0をシフトさせたときの,推定値からの残差

\(\Large \displaystyle Se = \sum_{i=1}^{n} \left[ y_i -\hat{a_0} \ Exp(- \hat{a_1} \ x_i) -a_2 \right]^2 \)

であり,a1を,δ,だけシフトさせて,固定し,その際のa0, a2の推定値をソルバーで推定しました.

dS,を見ていただけるとわかるように,推定値,Seが一番小さく,今回は,ほぼ左右対称に増加していることがわかります.

グラフ化すると,

のように,二乗+定数できれいに近似できます.

二次曲線の近似においてもきれいに近似でき,

\(\Large \displaystyle y = 0.2531 + 1.2105 \ \delta^2 \)

ここで,分散値は,

・分散

\(\Large \displaystyle Ve = \frac{1}{n-3} \sum_{i=1}^{n} \left[ y_i -\hat{a_0} \ Exp(- \hat{a_1} \ x_i) \right]^2 = \frac{Se}{n-3} = \frac{0.2531}{6-3} = 0.08438 \)

であり(a0,a1,a2の二つのパラメータが3つあるので,自由度は,n-3),

\(\Large \displaystyle 1.2105 \ \delta^2 = 0.08438 \)

となるδがSEとなるので,

\(\Large \displaystyle \delta^2 = \frac{ 0.08438}{1.2105} =0.06971 \)

\(\Large \displaystyle SE_{a_0} = \sqrt{\delta^2} =\color{red}{0.2640} \)

と推定できます.

 

・Rによる推定

Rでの近似を行ってみると,

プログラムは,

xx <- c(0,2,3,4,6,9)
yy <- c(11,5,3,2,1.5,1.1)
plot(xx,yy)
fm<-nls(yy~a0*exp(-a1*xx)+a2,start=c(a0=10,a1=0.5,a2=1),trace=TRUE)
summary(fm)

で,結果は,

Parameters:
Estimate Std. Error t value Pr(>|t|)
a0 10.16402 0.37708 26.954 0.000112 ***
a1 0.49456 0.04408 11.219 0.001518 **
a2 0.88393 0.26931 3.282 0.046348 *

となり,Kyplotにおいても,

推定値 標準誤差(SE)
A1 10.16401 0.377081
A2 0.494562 0.044083
A3 0.883926 0.269308

と同じ結果となり,今回の推定値と,ほぼ一致,します.

微妙に異なるのが気になりますが....

「統計解析の初歩」,の「1.2 非線形最小2乗法の基本的な考え方」には,

δのずらした値0.5 を変えると結果は異なり、近似標準誤差の精度が変わる
統計パッケージにより、近似標準誤差の値は幾分異なる

とあります.ここで,”0.5”,がどこから出てきたかはわかりません.そもそも横軸(x軸)の範囲に依存しちゃいますし..

そこで,σをとる値(±での平均値)でどう推定値が変わるかを調べてみました.その結果が,

とどのδでもRなどの推定値より下回っていました....謎です...

次ページに今までのまとめを記します.

 

l tr